Long range outward migration of giant planets, 
with application to Fomalhaut b 
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ABSTRACT 

Recent observations of exoplanets by direct imaging, reveal that giant planets orbit at a few dozens to 
more than a hundred of AU from their central star. The question of the origin of these planets challenges 
the standard theories of planet formation. We propose a new way of obtaining such far planets, by outward 
migration of a pair of planets formed in the 10 AU region. Two giant planets in mean motion resonance in 
a common gap in the protoplanetary disk migrate outwards, if the inner one is significantly more massive 
than the outer one. Using hydrodynamical simulations, we show that their semi major axes can increase by 
almost one order of magnitude. In a flared disk, the pair of planets should reach an asymptotic radius. This 
mechanism could account for the presence of Fomalhaut b ; then, a second, more massive planet, should be 
orbiting Fomalhaut at about 75 AU. 

Subject headings: planetary systems: formation — planetary systems: protoplanetary disks — methods: numerical 



Introduction 



Most of the known exoplanets are gaseous giants, with 
semi-major-axes below 2 AU, or even as close to their 
parent star as 0.01 AU. Inward planetary migration is 
generally considered as a necessary phenomenon to ex- 
plain this distribution, because these plan ets should not 
form there. In the core-accretion model ( Pollack et al.1 
19961 ) , giant planets should form beyond the ice-line (lo- 



star (|Kalas et al.ll2008t iMarois et al.ll2008l ). It has been 
proposed that these large semi-major-axes may be due to 
scattering with other giant planets after formation in the 
5-20 AU region (e g. lScharf fc MenovJl2009l: IVeras et al " 



20091) . iBolevi (l2009h . IClarkd (|2009h . and lRafikovl (l2009h 



cated at about 4 AU for a solar type star). The general 
expectation is that planets form at distances comparable 
to those characterizing the orbits of the giant planets of 
our solar system, because further out the dynamical (and 
accretional) timescales are too long. 

However, giant planets have been recently observed by 
direct imaging at a few dozens to 120 AU from their host 



sugg est that planet forma tion by gravitational instability 
(see [5 urisen et al.l 120071 for a review) can be effective 
beyond ~ 50 AU, because the cooling time relative to 
the dynamical time becomes short, and the Toomre Q 
parameter may be smaller in the colder, outer parts of 
the disk. 

In this Letter, we explain how resonant interactions 
between two planets in a common gap can lead to out- 
ward migration, possibly explaining the presence of these 
planets. We outline the mechanism in Sect. [51 and 
present proof-of-concept simulations in Sect. [3] The ef- 
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fects of disk structure are discussed in Sect. [3J where we 
show that a pair of planets can be driven to an equi- 
librium radius in a flared disk. Caveats are discussed 
in Sect. and our main conclusions and applications of 
this mechanism, paying special attention to the case of 
Fomalhaut, are presented in Sect. [5] 

2. Outward migration in resonance 

Many authors have observed that two planets migrat- 



resonance (e.g. Snellgrove et al.1 2001 


: Papaloizou 2003: 


Kiev et al.l 2004 Pierens & Nelson 


2008; Crida et al. 


2008). iMasset & Snelkrovd (2001) found that if the 



outer planet is significantly less massive than the inner 
one, and if the two planets open overlapping gaps, then 
the migration of the pair can proceed outwards. In this 
situation, the inner planet lies close to the inner edge 
of the common gap, and far from its outer edge ; there- 
fore, it feels a positive torque from the inner disk, and 
no torque from the outer disk. Symmetrically, the outer 
planet mostly feels a negative torque from the outer 
disk. If the inner planet is more massive, it feels a larger 
torque, in absolute value, than the outer planet, and the 
total torque applied to the pair of planets is positive. 
The pair of planet therefore moves outwards. 

In order for this process to last on the long term, the 
material lying outside of the common gap must be fun- 
neled towards the inner disk. Otherwise, it piles up at 
the outer edge of the gap, which eventually reverses the 
torque balance. In addition, the inner disk must be re- 
filled. In fact, assuming a local damping of the outer 
planet's wake, one finds in the simulations an outwards 
drift rate of the disk material significantly smaller than 
that of the planets. Thus, the gas just outside the gap 
is necessarily caught up with the outer separatrix of the 
outer planet, and thrown inwards. 

The whole process may be enhanced by the horseshoe 
drag due to gas of the outer disk pas sing through the 
comm on gap - like in type III migration (jMasset fc Papaloizoul 
20031 ). However, the outwards migration of the pair is 
observed at all disk mass in numerical simulations. Since 
the corotational effects like in type III migration become 
virtually negligible at small disk mass, it is primarily the 
wake torque imbalance that drives the outward migra- 
tion. 

Another clue of the importance of the Lindblad 
torque imbalance on the process is the sensitivity of 
the drift rate on the disk's aspect ratio H/r (H be- 
ing the scale height of the disk, and r the distance 
to the star), as the o ne-sided Lindblad torque is pro- 



portional to (H/r)- 3 . IMasset fc Snellgrovel (|200lh and 



This effect is further enhanced by the fact that the tidal 
truncation of the outer edge of the common gap by the 
outer planet is less pronounced in thicker disks. 



3.1 



Long term simulation 
Code presentation 



Accurate numerical simulation of the evolution of 
giant planets requires simultaneous computation of 
planets-disk interactions and global evolution of t he disk. 
Therefore, we use the 2D1D version of FARGO (|Massed 
fcOOOal lbM in which the standard 2D polar grid is sur- 
rounded by a ID grid to compute the disk evolution on 
all its physical extension ( Crida et al. 20071 ). 

The 2D grid extends from 0.45 to 12.61 L, where L is 
the length unit (taken as 10 AU). The ID grid - that is, 
the disk- extends from 0.05 to 20 L. The resolution is 
5r/r = 0.01, constant in the two grids. In the 2D grid, 
59 = 5r/r. 

3.2. Disk and planets settings 

The inner regions of protoplanetary disks are poorly 
constrained, and disk profiles must be assumed for typi- 
cal inner disk migr ation studies. In c ontrast, millimeter 
interferometry (e.g. iPietu et al.l l2007) places strong con- 
straints on disk structures between 30 and 200 AU. 

Of particular relevance to scenarios of migration is 
the temperature profile, which can be recast as H/r = 
c s /rf2, where c s is the sound speed and Q the angular 
velocity. The aspect ratio follows H/r = ho x (r/Vf . 
The flaring ind ex (3 is typic ally 0. 25 and range s from 0.17 



to 0.32 (jPietu et al.ll2007t ). The IPietu et al.1 results are 



compatible with ho = 0.045, for fj, = 2.4 (mean molecular 
weight) and T = 30 K at 100 AU. Because we use a 
locally isothermal equation of state in our simulations, 
the temperature and aspect ratio are a function of radius 
only, independent of time. 

For the surface density profile, we assume E(r) = Eo x 
(r/L)~ p . Detection of dust in the continuum constrains 
the profile of the dust column density multiplied by its 
emissivity (hence its temperature). Assuming a uniform 
dust to gas ratio, this yields 1 — 20 + p m 2, thus p m 1.5 , 



value adopted in this work (see e.g. IPietu et al.ll2007h . 
Eo is much less constrained. However, unless the disk 
mass is large enough for type III migration to affect the 
dynamics, our results should be essentially independent 
of So. More precisely the drift rate should scale with 
Eo _1 ; this is approximately what is seen in our numerical 
simulations with Eo varying by one order of magnitude. 
We therefore adopt values of the surface density that offer 
a good trade off between speed and the onset of type III 



Morbidelli k. Cridal (|2007l ) find that the migration rate 
is a decreasing function of the aspect ratio of the disk. 



J http : //f argo . in2p3 . f r 



2 



migration of the outer planet that would take it away 
from the resonance with the inner one. 

Calculations are run with a IShakura fc Sunvaev 



(1973) a-coefficient between 0.001 and 0.01 in order 
to explore the importance of the turbulent viscosity 
v = ac s H . 

The two planets have masses Mi = 3 x 10~ 3 M* 
and M 2 = 10~ 3 M*, i.e. 2M Jupi ter if Af* = 2M Q . 
The planets are initially on circular orbits at a\ = L, 
ci2 = 2L, and do not feel the disk potential for the 
first 100 (Sect. and BT3J or 400 (Sect, gj) orbits of 
the inner planet, so the disk can adapt to their pres- 
ence. Denoting r H = r(M p /3M^) 1 ^ 3 the Hill radius of a 
planet, the region inside 0.6 rn for a given planet is not 
taken into account when calculating the force of the disk 
on the plan e t, usi ng a smooth filter, in agreement with 



Crida et all (|2009h . 



3.3. Result 

Consider the case a = 0.01, H/r = 0.045(r/L) 1 / 4 , 
and E = 1.5 x 10~ 3 M„L" 2 . Assuming M* = 2M Q 
(like Fomalhaut) and L = 10 AU, this is 6 .6 times 
dense r than the Minimum Mass Solar Nebula ( HavashJ 
Il98l[ ) ; the disk parameters are summarized in Table [TJ 
The migration path of the planets in this disk is plot- 
ted in the top panel of Fig. [TJ After the release, the 
outer planet migrates rapidly inwards, and is captured 
in 2:1 mean motion resonance. Then, the two planets 
migrate smoothly outwards together. The resonant an- 
gle <7i = — Ai +2A2 — 0J1 is displayed in the bottom panel. 
The outer planet reaches 10 L (100 AU) in 4 x 10 5 years, 
and the migration speed decreases as the planets move 
further from the star. At distances larger than 100 AU 
for the outer planet, the orbits seem to converge to a 
stable equilibrium radius. 

Assuming that the migration speed is proportional to 
the gas density, the outer planet would reach 100 AU in 
4 million years in a disk that is only two thirds the mass 
of an MMSN. Because this time is roughly the lifetime 
of a protoplanetary disk, a disk that is 10 times lighter 
than considered above represents the minimum mass re- 
quired to migrate the outer planet beyond 100 AU. Such 
a disk would be light for a 2 M© star, so we conclude 
that there is sufficient time for the resonance mechanism 
to transport a massive planet to more than 100 AU. 



Table 1: Parameters of the disk used in Sect. 13.3 



r 


T [K] 


H/r 


S [kg.m~ 2 ] 


Q 


10 AU 


104 


0.045 


2666 


9.55 


20 AU 


73.5 


0.0535 


943 


8.0 


50 AU 


46.5 


0.0673 


238 


6.4 


100 AU 


32.9 


0.08 


84.3 


5.37 



4. Extension of concept 

To understand this process in general, we have per- 
formed additional short-term simulations, using a smaller 
2D grid (extending from 0.4 to 5 L), with [3 = 0.25, 
p = 1.5, S = 1Q- 3 M*L- 2 (1778 kg.m~ 2 in the case 
M* = 2M , L = 10 AU), and with varying ho and a. 
In each case, the migration speed of the pair of plan- 
ets is measured after the resonance capture. The result 
is displayed in Fig. [2] The figure is almost filled with 
symbols, which means that outward migration is possi- 
ble for a wide range of parameters. More precisely, the 
starred symbols (connected with the vertical solid line) 
display the migration rates obtained with aspect ratio 
H/r = 0.06{r/L)^ 4 , for various values of a between 
10~ 2 (bottom) and 10 -3 (top) ; dln(a)/df is a decreasing 
function of a, and there exists a critical a c such that the 
migration is directed outward for a < a c and inward for 
a > a c . 

Reciprocally, for a fixed a, dln(a)/d< is a decreasing 
function of H/r (see the short-dashed (a = 10 -2 ) and 
long-dashed (a = 10 -3 ) curves). In both cases, there ex- 
ists a critical aspect ratio h c such that the migration is 
directed outwards for ho < h c and inwards for hg > h c . 
This explains partly the slowing down of the migration 
as the semi major axes increase in Fig. [TJ and suggests 
that in a flared disk with uniform a, the pair of planets 
should reach the location where H/r = h c . This is a sta- 
ble equilibrium point as far as migration is concerned : 
in the inside, H/r < h c , and migration is directed out- 
wards ; further from the star, H/r > h c and migration is 
directed inwards. If H/r increases with r, there should 
be convergence toward the place where H/r — h c . The 
evolution of the density profile also plays a role, but the 
study of this parameter is beyond the scope of this letter. 

On e can expect that the efficiency of the lMasset fc Snellgrovei 
( 200l[ ) mechanism is directly related to the size of the gap 
opened by the outer planet. To demonstrate this, we plot 
in Fig.[3]the same points and curves as in Fig. [2l but the 
migration speed is expressed in viscous time t v — r 2 /v, 
and the x-axis reports the parameter P, defined as : 



3_ff 

4 r H 



(1) 



This parameter is precisely related to the width an d 
depth of the gap of the outer planet (jCrida et al. I l2006l) . 
The strong correlation is obvious. For a broad range of 
H/r and a, the dependence of the migration speed on 
these two parameters can be approximated by a depen- 
dence on the sole parameter P. The migration rate is 
close to stationary for 2 < P < 2.5, whatever the values 
of a and H/r that combine into this value of P. 

Then, using Eq. |TJ), P = 2.5 gives the critical aspect 
ratio that halts the migration as a function of a, or a c for 
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5x10 10 1.5x10 2x10° 2.5x10° 3x10 3.5x10 4x10 
time [years] 

Fig. 1. — Top: Semi major axes a\ and 02 as a func- 
tion of time (assuming M* = 2M©). Bottom: Resonant 
angle <j\ = —\\ + 2X2 — u>i as a function of time. 
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Fig. 2. — Migration speed of a resonant pair of planets, 
for various values of a and ho- + symbols: a = 10~ 2 . 
x symbols: a = 1CP 3 . Starred symbols: ho — 0.06, and 
a = 10~ 3 , 2 x 10~ 3 , 5 x 10~ 3 , 10~ 2 from top to bottom. 



a given H/r. The top panel of Fig. |4] shows the value of 
the aspect ratio that makes P = 2.5 for the outer planet, 
as a function of a, for 2 x 10~ 4 < a < 0.2. The radius 
at which the outer planet would stop is then 



(H/r) 



1//3 



L 



This is displayed in the bottom panel for ho = 0.06, so 
that L can be understood as the radius in the disk where 
H/r = 0.06. The equilibrium radius r c depends on a, j3, 
and Mi- Since (3 is generally small, r c can be very large 
with respect to L. 

Note that for a given value of P, the migration speed 
may depend on other parameters. Therefore, the critical 
value for stationary migration P w 2.5 holds only in the 
case considered here, that is with M\/Mz — 3, (3 — 0.25, 
S = 10~ 3 Af*L~ 2 and p = 1.5. A full analysis of all the 
parameter space is beyond the scope of this letter. 

5. Discussion 

5.1. Gas accretion on the planets 

A shortcoming of the mechanism presented here is 
that it disregards the accretion of gas onto the plan- 
ets. As the pair of planets proceeds outwards, the gas 
proceeds inwards through the common gap. The outer 
planet should presumably accrete most of the incoming 
material, narrowing the mass difference with the inner 
planet, so that the torque balance could cancel or re- 
verse at some point. Accretion of the nebular gas onto 
a giant planet is a complex process. It is dependent on 
the global structure of the flow, on the local properties of 
the circum-planetary disk or envelope (poorly described 
in numerical simulations), and on the micro-physics of 
the nebula, such as the opacity of the grains. Most de- 
tailed studies of gas accretion onto giant planets can't 
be applied to the present case, because the conditions 
are different from the conditions that prevail at 10 AU, 
where these studies apply. This appeals for a thorough 
investigation of gas accretion in the conditions of outward 
migration in the outskirts of the disk, which is beyond 
the scope of this work. 

5.2. Self-gravity of the gas and 3D effects 

The gas self-gravity may yield, depending on the equa- 
tion of state of the fluid and the mass of the disk, to 
a vertical com pression of the fluid i n the vicinity of 
the shock (e.g. Wev fc Durisenl [20061 ) . This would in- 
crease the Lindblad torque, more efficiently on the side 
of the most massive planet, thus it would actually en- 
hance the mechanism of outward migration. But 3D 
effects apply mostly at high order resonances, close to 
the planet, which are depleted in our case. In addition, 
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2 3 
P = (3/4)H/r H + 50/qR 



Fig. 3. — Migration speed of a resonant pair of planets, 
normalized to the viscous time, as a function of P. Same 
key as in Fig. [2] 
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thr ee dimensiona l semi -analytic calculations performed 
bv lTanakaet all (l2002h show that the angular momen- 
tum is essentially exchanged with waves that have no 
vertical structure. Therefore, owing to the size of the 
gap that our pair of planets carves in the disk, and ow- 
ing to the relatively large value of the Toomre parameter 
(see Table [T]), we feel confident that non-self-gravitating, 
two dimensional isothermal simulations capture the main 
features of our mechanism and provide an acceptable or- 
der of magnitude of the outward drift rate. 

5.3. Equation of state 

In the above simulations, the equation of state was lo- 
cally isothermal. This is realistic at large distances from 
the star, where the opacity is low and the cooling time is 
short with respect to the dynamical time. However, at 
~ 10 AU from the star, the cooling time can be larger 
than the orbital period. This can affect the aspect ratio, 
because of the heating of the disk by the wakes launched 
by the planets. 

Additional simulations have been run, with the set- 
tings of Sect. I3.3[ an d computing the full energy equa- 



tion, with viscous heating Q+, and no heat diffusion 
in the disk plane, but vertical radiative cooling Q_ = 
2ctrT 4 / kT, (where an is the Stefan-Boltzmann constant). 
The opacity n is constant and tuned for the unper turbed 
disk to be in thermal equilibrium (see Eq. (6) of Cridal 
2009). This gives a cooling time of 111 fl 1 initially. The 
aspect ratio increases to about 0.06 in the planets region 
after 100 orbits of the inner planet. This makes the mi- 
gration rate after the resonance capture to decrease to 
~ 10~ 4 AU/year instead of ~ 4 x 10" 4 AU/year. 

With a realistic opacity, given by iBell fc Linl (|l994 
Appendix), T and H/r are increased for r < 10 AU, 
and shrink for r > 15 AU because we don't take heat- 
ing from the star into account. At 10 AU, the cooling 
time is then ~ 400 f2 _1 . Migration proceeds outwards 
at a speed initially similar to the locally isothermal case 
(which was expected as Hj r is not much changed around 
10 — 15 AU), slightly accelerating while reaching the outer 
regions where H/r is smaller. 

So, it seems that migration can proceed outwards on 
the long range also in non locally isothermal disks, even 
if a self-consistent radiative disk model still has to be 
tested. 



Fig. 4. — Top : Solutions of the equation P = 2.5 for 
three values of M2/M*. Bottom: Corresponding r c , for 
h = 0.06. 



6. Applications 

We find that a pair of planets formed in the 5 — 20 AU 
region can reach large semi major axes, if the outer planet 
is lighter than the inner one and they orbit in resonance 
in a common gap in the protoplanetary disk. The pair 
of planets tends to reach regions in the disk where H/r 
has a critical value (of the order of 0.1 — 0.15, depending 
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on the other parameters) ; this can be ten times further 
from the star than where they formed. After the disk has 
dissipated, the evolution of the planets may lead to the 
disruption of their original resonant configuration, due 
to various processes : (i) the eccentricities can rise too 
much if the planets are too massive ; (ii) a third planet 
could destabilize the system at some point ; (iii) the outer 
planet may enter a debris disk, in which it migrates in- 
dependently of the inner planet by planetesimals scatter- 
ing. Such a phenomenon may account for the HR8799 
system, where the outermost planet is the lightest one 



but no resonance has been identified so far (Marois et al 



2008) 



Otherwise, the planets should remain on orbits in reso- 
nance far from the star. Their eccentricities should there- 
fore be non zero but moderate : during the outward mi- 
gration, e is excited by th e resonance and si multaneously 
damped by the disk (e.g. ICrida et al. 20081 ). In the sim- 
ulation presented in Fig. [TJ the eccentricities are about 
0.01 and 0.02 — 0.03 for the inner and outer planet respec- 
tively. In contrast, the scattering of giant planets leads 
to very eccentric orbits (up to e = 0.8). The typical final 
eccentricity of planets formed by gravitational instability 
is not known. 

The mechanism presented here could account for the 
case of Fomalhaut b. Indeed, this p lanet has a mass 
smaller than 15 x 1 Q- 3 M„ ( 3Mj UDi t r ) (iKalas et al.l l2008l: 
Chiang et all l2009h . while iQuillenl ((2006J) predicted a 
planet smaller than Saturn. Thu s, a more ma s sive, i nner 
planet could exist. In addition, IKalas et al.l (|2008l ) no- 
tice that Fomalhaut b is not apsidally aligned with the 
dust belt, which suggests the presence of additional per- 
turbers in the sy s tem. Observations in the M-band by 
Kenworthv et al. ( 20091 ) rule out the presence of a planet 
more massive than 2 Afj up i tC r in the range 8 — 40 AU from 
Fomalhaut. However, if our mechanism occurred in the 
Fomalhaut system, then the second planet should be or- 
biting at about 75 AU from the star, with a mass of the 
order of 1 to 10 Jupiter masses. 

Simulations have been performed on the hyades cluster of 
the D.A.M.T.P., on the on the hpc-bw cluster of the Rechen- 
zentrum of the University of Tubingen, and on a 92 core clus- 
ter funded by the program "Origine des Planetes et de la Vie" 
of the French Institut National des Sciences de I'Univers. A. 
Crida acknowledges STFC and J. Papaloizou. We acknowl- 
edge enlightening discussions with V. Pietu. We also thank 
the anonymous referee for his/her insightful remarks, which 
improved the paper. 
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